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Monotone Solutions of a Nonautonomous 
Differential Equation for a Sedimenting Sphere 

CNl" 

q ' Andrew Belmonte, Jon Jacobsen, & Anandhan Jayaraman 

Abstract 

We study a class of integrodifferential equations and related ordinary 
04 ' differential equations for the initial value problem of a rigid sphere falling 

through an infinite fluid medium. We prove that for creeping Newtonian 
flow, the motion of the sphere is monotone in its approach to the steady 
state solution given by the Stokes drag. We discuss this property in terms 
of a general nonautonomous second order differential equation, focusing 
on a decaying nonautonomous term motivated by the sedimenting sphere 
problem. 
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1 Introduction 

A rigid sphere falling through a viscous medium is a classic problem in fluid 
dynamics, which was first solved in the steady state for the limit of vanishingly 
small Reynolds number in an infinite domain by G. G. Stokes in 1851 E3). 
The time-dependent approach to the steady state allows the partial differential 
equation for the sphere and fluid to be reduced to an integrodifferential equation 
for the motion of the sphere in an infinite medium. The physical effects included 
in this equation are the buoyancy which drives the motion, the inertia of the 
sphere, the viscous drag, an added mass term, and a memory or history integral 
0- 

In the case of a Newtonian fluid, the main effect of the memory integral on 
the dynamics is to modify the approach to steady state from exponential to 
algebraic. The integral also makes the equation effectively second order, though 
it is generally accepted that no oscillations occur as the sphere reaches its steady 
state value jj], ^5| . Physically it is clear that no oscillations can occur due to 
the absence of a restoring force against gravity, and a sphere released from rest 
in a Newtonian fluid at low Reynolds number is observed to reach its terminal 
velocity monotonically. However, it is not directly evident mathematically that 
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Figure 1: Motion of a 1/8" diameter teflon sphere falling through an aqueous 
solution of 6 mM CTAB/NaSal: a) position vs time; b) calculated velocity vs 
time. 



oscillating solutions are precluded, particularly as the governing integrodiffer- 
ential equation can be transformed to a nonautonomous second order ordinary 
differential equation which has the form of a harmonic oscillator j?], |24| . 

In a non-Newtonian fluid, such as a polymer solution, a falling sphere is 
often observed to undergo transient oscillations before reaching its terminal 
velocity . These oscillations occur due to the elasticity of the fluid, which 

provides a restoring force . The steady state value is of primary concern in 
many applications, and much work focuses only on this aspect of the problem. 
The oscillations which occur during the approach to steady state have been 
reproduced in a linear viscoelastic model by King and Waters jl5| . 

More recently, nontransient oscillations of falling spheres (and rising bub- 
bles) have been observed in specific aqueous solutions of surfactants (wormlikc 
micellar solutions) J3], These observations were initially made for a bubble 
in the wormlike micellar fluid CTAB/NaSal [2l[], which showed oscillations 
in its position and shape. The shape oscillations included an apparent cusp 
which momentarily appears at the trailing end of the bubble. Such a cusplike 
tail is a well known property of rising bubbles in non-Newtonian fluids JiT| , |l7| . 
which we initially believed to play an important role in the micellar oscillations. 

Subsequent observations of solid spheres which also oscillate while falling 
through the same solutions made it clear that the cusp is not involved in the 
phenomenon, and that another explanation must be sought. Unlike a sedi- 
menting sphere in a conventional non-Newtonian fluid, these oscillations do not 
appear to be transient ||l3f . An example is given in Figure 0, which shows the 
motion of a 1/8" teflon sphere falling through a tube (L = 98 cm, R — 3.2 cm) 
filled with a 6 mM 1:1 solution of CTAB/NaSal ||. 

Our attempts to model this phenomenon brought to our attention the un- 
usual aspects of the integrodifferential equation for a falling sphere. We prove 
here that the equation for sedimenting sphere in a Newtonian fluid in the limit of 
zero Reynolds number (creeping flow) does not admit oscillating solutions, de- 



EJDE-XXXX/XXAndrcw Bclmonte, Jon Jacobsen & Anandhan Jayaraman 3 



spite some appearances that it does. This result is due to the special properties 
of the error function when multiplied by oscillating functions. It is ultimately 
related to the stability of nonautonomous ordinary differential equations with 
monotone secular terms, which is appropriately viewed as an initial value prob- 
lem, and not in terms of linear stability analysis around the terminal velocity. 



2 The Motion of a Sedimenting Sphere 

We begin by reviewing some classical results for the equation of motion govern- 
ing a falling sphere in a viscous Newtonian fluid of infinite extent (no bound- 
aries) . 



2.1 Equation of Motion of the Sphere 

An incompressible fluid in the absence of body forces is described by the equa- 
tions 



' du ,_ 
P 1 or + ( M ' V ) u 



dt 



divw 



= — Vp + div er, 
- 0, 



(2.1) 
(2.2) 



where p(x, t) is the density of the fluid, p is the pressure, u(x, t) is the velocity 
field for the fluid, and <j(x, t) is the extra stress tensor, which measures force 
per unit area (other than pressure) in the present configuration of the fluid. A 
Newtonian fluid is a fluid for which the stress tensor a is linearly related to the 
rate of strain tensor D through the relation 



a = 2uD, 



(2.3) 



where u is the viscosity of the fluid and D = (Vu + (Vu) T )/2 is the symmetric 
part of the velocity gradient VS. From (JO), (2. 2), and (2.2) one obtains the 
Navier-Stokes equation: 



du 

at 



(u-V)u 



-Vp + fiAu. 



(2.4) 



Non-Newtonian fluids are fluids for which the assumption (2^3) is invalid. 
For instance, polymeric and viscoelastic fluids often fail to conform to the in- 
stantaneous relation between stress and velocity gradients implicit in (2.3). In 
general <j = a(t, D) will depend nonlinearly on D and on the past history of 
stress in the fluid. 



By choosing a time scale and an appropriate length scale, (2.4) can be written 
in a nondimensional form 



du 



Re (u ■ V)u = -Vp + Aft 



div ft 



0. 



(2.5) 
(2.6) 
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where, u and p are the nondimensional velocity and pressure. The dimensionlcss 
constant Re is called the Reynolds number and it measures the relative impor- 
tance of inertial effects to that of viscous effects. When the inertial effects are 



negligible (Re = 0), equation (2.4) is called the Stokes equation. In this paper 



we restrict our analysis to this situation. 



Stokes solved the steady version of (2. 5)- (2. 6) for the case of sphere falling 



though the fluid for vanishing Reynolds number. The Stokes solution gives the 
steady state drag on the sphere of radius R falling through a fluid with a steady 
speed U to be F — 6irpRU. However, in order to solve the transient problem of 
falling sphere, we first solve the problem of sphere oscillating with a frequency 
uj and compute the drag experienced by the sphere as a function of u>. The drag 
experienced by a sphere falling at a arbitrary speed U(t) can then be computed 
as a Fourier integral of this drag. 

Consider a sphere of radius R and density p s in a Newtonian fluid of den- 
sity p and viscosity p. The exterior Stokes flow driven by small oscillations 
of the sphere at a frequency u> can be solved exactly E6j, which leads to a 
hydrodynamic force dependent on both U and dU / dt : 

F = 6npR U+*)U + 3nR 2 p6 (l + ^f) % (2.7) 



6 J r V 96 J dt 



where 5 = y/2vjui is a diffusive lengthscale common to Stokes problems, and 
v = p/p is the kinematic viscosity. Using this, the general time-dependent 
problem of the motion of a falling sphere can be reduced from a partial to an 
ordinary differential equation for the speed U(t) of the sphere, an exact equation 
which takes into account the motion of the surrounding fluid j^, ||| . For a sphere 
moving with an arbitrary speed U(t), the hydrodynamic drag it experiences can 
be calculated by representing U(t) as a Fourier integral: 

y* oo 

U(t) = / U u e- iwt duj. 



The drag for each Fourier component is then given by (2.7). The total hydrody- 
namic drag on the sphere is obtained by integrating over all Fourier components, 
leading to 

F drag = 6irpRU(t) + \ P V^- + 6irpR\[^ f -^L ds (2.8) 
2 dt V t J-oo yt - s 

where the first term represents the steady state drag on a sphere falling with a 
velocity U, the second term represents the added mass term (the force required 
to accelerate the surrounding fluid) , the third term is the Basset memory term, 
and V = 4ttR 3 /3 is the volume of the sphere. If the sphere starts from rest, then 
the lower limit of the integral in ( |2.8| ) starts from r = instead of r = — oo. The 
expression for the unsteady drag force can then be substituted into the balance 
of force equation for the sphere: 

PsVU (t) Ffj UO y Fd ra g. 
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Thus the equation of motion for the sphere is 

(Ps + ^)VU'{t) + QnpRU{t) + 6ttR\[^ [ -^2fL ds 
2 V T Jo V*-« 

= (Ps - p)Vg, (2.9) 
which can be rewritten in the simpler form 

U'(t)+BU(t)+Q [ -^M=ds = M, (2.10) 
Jo Vt-s 

where 

B = wwtnry < 2 -"» 

« " W^S' (2 ' 12 > 

M = (2.13) 

and Ap = /3 S — p is the density difference which drives the motion. In this 
approach the motion of the sphere is described by an integrodiffcrential equation 
whose integral term has the same singularity as Abel's equation Jli], 01. Note 
that this equation is only valid in the limit of zero Reynolds number [ fl9| , |l8| . 

Physically one expects the solution U(t) to approach a terminal velocity. It 
is clear from (2.10) that the only steady state solution (U' = 0) possible is 

which is the classical result of balancing the Stokes drag with the buoyancy. 

2.2 Solving the Integrodifferential Equation 

We first rewrite the integrodifferential equation (IDE) for the sphere in a nondi- 
mensional form using Uo as the velocity scale, and 1/B, the viscous diffusion 
time, as the time scale. The variables are 

r = Bt and u(t) = U(t/B)/U . 



With this rescaling (2.10) becomes 



u'(t)+u(t) + *[^ f-^=da=l, (2.15) 
V * Jo Vt-s 



where the control parameter n is given by 



ttQ 2 _ 9p 
~B~ ~ 2p s +p 



(2.16) 
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Thus the motion of the sphere depends only on the relative densities of the 
sphere and the fluid through the parameter k. The density of the sphere p s can 
range from zero to infinity, which implies a parameter range < k < 9. When 
the density of the sphere is equal to the density of the fluid, k = 3 and Uq = 0. 
Here we will only be concerned with falling spheres, for which < n < 3. 



Although we will solve the IDE (2.15) directly, it is of interest to connect the 



problem to ordinary differential equations (ODE's) and discuss some important 
consequences therein, especially with regard to the st abilit y of the terminal 
velocity solution. Following Villat p4j , we can rewrite fl2.15 ) as an ODE using 
Abel's Theorem (see Appendix |a|) : 

u" + (2 - k)u' + u = 1 - ,/— , (2.17) 

V 7TT 

tt(0)=0, «'(0) = 1. (2.18) 

More general initial conditions U(0) j= and U'(0) = M - BU(0) lead to a 
slightly different ODE: 



u" + (2-k)u' + u = 1 + J — (u(0)-l), 



u(0) 



t u'(0) = l-t 



(2.19) 
(2.20) 



Note that it follows from equation (2.15) that u'(0) is prescribed in terms of 
u(0), and thus the second order ODE we have obtained requires only one initial 
condition. 

Since we are investigating the possibility of steady-state oscillations of a 
sedimenting sphere, we are primarily concerned with the asymptotic behavior 
of (2.17). Moreover, since the nonautonomous term tends to zero as t —> oo, one 
might expect the stability of (2.17) to mimic the homogeneous problem. With 
this in mind, let a and denote the roots of the characteristic equation 



rn 



(2-«)m + l = 0. (2.21) 

It is readily verified that ( 2.2l| ) has complex roots for < K < 4. Moreover, the 
roots have positive real parts for 2 < n < 4. Since the relevant range of k for a 
falling sphere is < k < 3, one sees that oscillations are not a priori precluded. 
In terms of the actual densities of the fluid and sphere the condition for complex 
roots corresponds to p s > (5/8)/?, which is true in the case of a heavy sphere 
falling through a lighter liquid (p s > p). If additionall y p s < (7/4)/), then the 
complex roots have positive real parts. If we rewrite ( |2.17 ) in the asymptotic 
limit (t — > oo) as a first order linear system 



x = y 

V = (l-x) + (K-2)y, 



(2.22) 
(2.23) 



where x = u, and y — u', then (x,y) = (1,0) is the unique equilibrium point, 
which corresponds to the terminal velocity. The eigenvalues of this system are 
precisely a and /3, whence the equilibrium point becomes unstable. 
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Nonetheless, as we will show, even in this range (2 < k < 3), the solution to 
the full equation (2.17) monotonically approaches the value 1, corresponding to 
the monotonic approach to the steady Stokes value (2.14) for the actual velocity 
U. Clearly the nonautonomous term continues to play a dominant role in the 
stability of (2.17), despite its algebraic approach to zero. 

To solve for u = u(t), we return to the IDE ( 2.15 ) and apply the Laplace 
transform in the case u(0) = 0: 



£{u}(s) = 



s(s + + 1) 



or £{u'}(s) = s£{u} 



s + v^v^s + 1 



Since a/3 = 1 and \fu + \f]3 = V^, we may express this last equation in the 
form 



£{u'}(s) 



a- 13 



VP 



Moreover, the identity 



£{e at ErfcvW} 



1 



implies 



u'(r) 



a- (3 



VoIe QT Erfcv / aT - V^e^Erfc^//^ 



Finally, since 



y/a / e at EricVaidt = 2\ - + ^e Qt Erfc\/oi + C, 

J V 7T V" 



we find the solution 

u(t) = 1 + 



a- 13 



""ErfcVoT e^ r Erfc v / ^r 



VP 



(2.24) 



(2.25) 



This solution to the IDE (|2JJ) is also the solution to the ODE Q2.17| )-( |2.1S| ). 
Applying transform methods to the more general set of equations defined by 
fl2.19D -(f2~20l) one finds the solution 



u(t) = (1 - e)uo(r) + e 



u(0), u'(0) = 1 - e, 



(2.26) 



where Uq(t) is the solution defined by (2.25). Note that the solution for arbitrary 
initial velocity e is a simple rescaling of the solution for the sphere initially at 
rest. It is not obvious from the form of u in (2.25) that the values approach 1 
monotonically. Let us first investigate the asymptotic behavior of this solution. 
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2.3 Asymptotic Approach to the Steady Stokes Solution 

Finding the asymptotic behavior of the solution uq is straightforward. We em- 
ploy the asymptotic expansion of the error function H 



e z Erfcz 



1 



1 • 3- • -(2771-1) 

(2z 2 ) m 



(2.27) 



as z — > oo, provided | arg(z)| < 37r/4. Since < arg(Vai) < 7r/2 for < K < 3, 
we see that asymptotically 



Erfc 



CVT 



Thus the product of the exponential term and the error function approaches 
zero in the limit t — > oo. Using this expansion in ( p. 25 ) we obtain 



lim uq(t) = 1 



lim U(t) = U . 



(2.28) 



As the solution for any initial condition is a rescaling of uq, we see this limit 
applies for all values of w(0). 

Although the asymptotics of this solution are clear, the transient solution 
has some unusual properties. Numerical simulation of the IDE (2.15), or even 
attempts to plot the analytic solution (2.25), eventually blow up at large r when 
k is in the unstable range. Clearly the cancellation between the exponentially 
growing and decreasing terms are quite sensitive to numerical errors. This is 
an indication that the product e aT Erfc-y/ar should be considered as a special 
function with its own properties. 



3 Monotonicity of the Transient Solution 

We begin with the transient solution to the IDE or ODE considered in the 
previous section. In addition to the insensitivity of the transient solution to the 
real part of the homogeneous roots, it is surprising that the nonzero complex 
part does not lead to any oscillations in the velocity of the sphere, although 
there has sometimes been confusion on this point regarding transient oscillations 



24 . Experimentally the sphere in a Newtonian fluid has never been observed 



to oscillate, in contradistinction to most non-Newtonian (p articu larly elastic) 



fluids [|l| |25|]. We will show that the solution u defined by ( 2.25 ) is monotone 
as a function of r for all k £ (0, 4). Although this may be well known, we have 
not yet found a reference to any proof other than the case a,/Jel f2q| , which 
corresponds to k > 4. 

Let us define the function Vi : C — > C by 



Vi(z) = e z ErfrVI = — / e~ s ds. (3.1) 



/7T 
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We shall refer to Vi as the Villat function, since this combination appeared in 
the explicit solution of the differential equation for the falling sphere problem 
by Villat (24j. Closely related to Vi(z) is the "plasma dispersion function" w, 
defined by 

2 

w(z) — e~ z Erfc(— iz). 

In fact, Vi(ai) = w{i^fod). Using the Villat function we may now prove the 
main theorem. 

Theorem 3.1 (Monotonicity) For each k S (0,4), the function 

' e aT Erfc^/cPF e 0r Erfcy/Jh' 



u{t) = 1 + 

a — p 

approaches the limit 1 monotonically. 



VP 



(3.2) 



Proof. We have shown that u(t) — > 1, thus it remains to show it does so 
monotonically. We will demonstrate this by proving u'(t) > for all t > 0. To 
this end, fix t > 0, k E (0,4) and recall that a and (3 denote the conjugate pair 
of roots of the polynomial m 2 + (2 — n)m + 1. Recall from ( 2.24 ) that 



u'(t) 



a- (3 



V^Vi(at) - V/3Vi(/3t) 



Since Erfc(z) = Erfcz, it follows that Vi(/3i) = Vi(ai) and 

KV^Vi(ai)} 



, M r- \fa Vi(ai) - y/a Vi(crf) / -' 
u (t) = \/k — = v K ' 



(3.3) 



Since 3{a} > for each k € (0, 4), it is evident from ( |3~3| ) that the sign of u'(t) 
is determined by the imaginary part of the function ^fa Vi(ai). For the plasma 
dispersion function w introduced above, the real and imaginary parts are given 
by (see e.g., § 7.4.13-7.4.14]) 



and 



3?(u>(x + iy)) = — 

7T 



^ /*°° 

f(w(x + iy)) = — 

7T 



ye 



(x — s) 2 + y 2 

(x — s)e~ s2 
(x — s) 2 + y 2 



ds (x€R,y>0) 



ds (xeR,y>0). 



It is readily verified that |a| = 1, thus in polar form we have a = e for some 



fixed 9 € (0, 7r), in which case Vi(crf) = w{iy/ai) = w(x + iy), where 



-Vt sir 



and 



Vi<, 
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Using this information we compute 

^s{^/a Vi(at)} = cos ^{w(x + iy)} + sin T^J $t{w(x + iy)} 

1 fO\ f°° ye-* 2 , . 

^■(ijLFW 4 (3 - 4) 

i 

7T \2J J_ oc (x - s) 2 +y 



Z C0S (o)/ 71 ^2 , „2 dS " ( 3 ' 5 ) 



In the last step we have used the fact that a;cos(6'/2) + ysin(0/2) = 0. Since 
i/a lies in the first quadrant, the prefactor of the last integral above is positive 
and we may conclude that ^{i/aV^at)} > provided 

se-" 2 f°° se~ s2 . . 

77; ds — / — 7. ds < 0. (3.6) 

,(x-s) 2 +y 2 ( s + Visinf) 2 +icos 2 f 1 7 

Let us denote the integrand as 



F(s) — — ~, r - where P(s) = ( s + y/ism - ) + t cos 2 -. 
P(s) V 2 / 2 

Note that P(s) > for s £ M. (recall t > is fixed). The proof is complete once 
the following two observations are made: 

(a) \F(-s)\ > F(s) for s > 0; 

W f R Fd8 = f™F(8)d8-J™\F(-8)\d8. 

To see (a), notice for s > we have < P(—s) < Pis), thus 

_ 2 _ 2 

P(-s) P{s) 

Observation (b) follows from a standard change of variables 

/> />0 />oo />oo />oo 

/ F{s)ds= / F(s)ds + / F{s)ds= / F(s)ds- / |F(-s)|ds. 

JK i-oo JO JO JO 



The two observations above imply the inequality ( |3.6| ) holds, in which case by 
( ft .3] ) and ( |3.5| ) we see it'(i) > 0. Since t > was arbitrary, the proof is complete. 





Corollary 3.2 The solution to the initial value problem (2.1i)-(2.2L) mono- 
tonically approaches its steady state value u = 1. 



EJDE-XXXX/XXAndrcw Bclmonte, Jon Jacobsen & Anandhan Jayaraman 11 



Proof. The proof follows from applying Theorem 3.1 to equation (2.26). ^> 



4 Related Aspects of the Newtonian Problem 

To investigate the generality of the above result, consider the nonautonomous 
linear damped harmonic oscillator equation for u = u(t) 

u" + bu' + u = l-G(t), (4.1) 

as an initial value problem with arbitrary initial conditions u(0) and u'(0). We 
are specifically interested in the case where G{t) — > as t —> oo, as opposed to 
the often studied case where G(t) is periodic (see e.g. {hJ). The Newtonian 
sphere problem ( [2.17 ) is a special case of ( |4.l| ), with b = 2—k and G(t) = ^/ ' n/nt. 



Making the change of variables v = u — 1 , we may simplify the equation to 

v" + bv'+v = -G{t) (4.2) 

so that v — solves the homogeneous equation. Note however that if G(t) ^ 0, 
then v = is not a solution to ( |4.2| ) for any t > 0. We are interested in 
the following question: what conditions on G{t) and b are necessary for the 
solution v(t) to remain monotone, even within the regime of instability for the 
homogeneous equation. 

As a first step in this direction we consider the following initial value problem 
for t > : 

v" + bv' + v = A (4.3) 

where b, A £ M. and to > are constants. The motivation for this form is to test 
the necessity of the singularity at t = in the monotonicity result of Section ||. 
To ensure complex roots, we assume b G (—2, 2) 



Using variation of parameters one finds a particular solution of (4.3) to be 
A 



v p(t) = -r, I V/? e " (t+to) (ErfiVa^ - ExIcJ a[t + t )) 

p — a L V / 

- Vae P{t+to} (Erfcv 7 ^ - Erfc^/3(i + i )) } , (4.4) 

where a and (3 are the roots of the characteristic polynomial m 2 + bin + 1. 
Employing the Villat function we may express equation ( |4.4| ) as 

v P {t) = ^-{ v /^Vi(at )e Qt -^Vi(/3t )e^}+^M(i + i ), 

where the function M is defined by 

M(t) = -1- {v^Vi(at) - V^Vi(/3<)} . (4.5) 
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In Section we proved M approaches monotonically for all b 6 (—2, 2). The 



general solution to (4.3) may be expressed as 



v(t) = C ie at + C 2 el n + v p {t) (4.6) 

Ci + ypAVi(at ) \ &at + r ^ _ Voavk^o) 1 e/3t 



+ AM(t + t ). (4.7) 

This equation clearly demonstrates how the long term dynamics of v(t) depend 
on the solution of the homogeneous problem. In particular, it shows that the 
solution will retain the stability properties of the homogeneous solution u nles s 



the coefficients C\ and C2 are chosen to zero out the first two terms in (4.7). 
The unique choice of C\ and C2 for this to happen are 

= _ Mf$_ and ^ = M/^ 

p — a p — a 



Moreover, it is clear from ( 4.7) that the solution in this case is v(t) — A M(t+to), 
with v(0) = AM(t ) and i>'(0) = AM'(t ). Thus, in this case, the solution is 
a translate of the monotone solution. The coefficients C\ and C2 are related to 
the initial conditions v(0) and v'(0) via 

Ci= m°w(q) and C2= mz^m. (4 .9) 

p — a p — a 

The values of v(0) = AM(tp) and v'(0) — AM'(t ) may also be obtained by 
solving equations ( |4.<j ) and ( [4.9|) . 

In summary, given b e (— 2,2), A > 0, and t > 0, for the equation 

A 

v" + bv' + v = (4.10) 

VMt + t Q ) 

there exists a unique choice of initial values v(0) = AM(to), v'(0) = AM' (to) 
such that the solution v(t) remains monotone for all t > 0. Therefore the 
presence of a singularity at t = in the nonhomogeneous term is not necessary 
to obtain a monotone solution. 

In light of the above analysis it becomes clear how the solution for the 
sedimenting sphere remains monotone in its approach to terminal velocity for 
all relevant values of k (i.e., sphere densities). From equations (2. 17 )-( p~18 ) 



we see that for each value of k, equation (4.10) describes the dynamics for 



the dimensionless velocity v = u — 1, with b = 2 — k, t = 0, and A = ^/k. 



Moreover, for each k £ (0,4) we have demonstrated that equation (4.10) with 
to = 0, b = 2 — k, and A = y/Ti, has a unique initial value for which the solution 
remains monotone, namely, 

,(0) = AM(0) = V^^^ = (4.11) 
a- p ^+ v/3 
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where a and (3 denote the roots of the polynomial m 2 + (2 — K)m + 1. However, 
since y/a lies in the first quadrant, \fu-V > 0, and the computation 

(Va + V/3) 2 = a + f3 + 2=-b+2 = n, 



together with ( 4.11 ), implies 

v(0) = -1. 

In other words, the particular relation between the parameters A and b decouples 
v(0) from all parameters, so that one obtains a monotone solution for all values 
of the sphere density. 

We conclude this section with a geometric interpretation of the monotonic- 



ity result. In particular, we focus on the interesting case of (4.3) when the 
parameter b E (—2,0). For these parameter values the solution v = of the 
homogeneous problem is unstable. We have shown that there is a unique set of 
initial conditions that defines a solution to the nonhomogeneous problem which 
approaches the unstable fixed point v — monotonically, despite the surround- 
ing instabilities. Thus we return to the fundamental puzzle posed in Section 



2.3: How is it that the nonautonomous term in (4.3), which decays to zero as 
t — > oo, can "stabilize" a trajectory for all t > 0, in the sense that this solu- 
tion approaches while all other trajectories diverge due to the instability of 
the linearized problem? The following observation resolves the puzzle. First, 
recall that the unique initial conditions for which the nonhomogeneous problem 
remains monotone are defined by 

v(0) = AM (to) and v'(G) = AM' (to). 

Second, note that as the amplitude A of the nonhomogeneous term tends to zero, 
the initial conditions (A M(i n ), A M'(t )) approach (0,0). This corresponds to 
the initial condition starting on the unstable equilibrium point, which is the 
unique initial condition for the homogeneous problem whose solution does not 
diverge. In other words we have a correspondence between the trajectories of the 
homogeneous equation and the nonhomogeneous equation, which is continuous 
with respect to the parameter A. The monotone solution is then the image of 
the unstable fixed point under this map. 



5 Summary and Conclusion 

In this paper we have studied the ODE model for a sphere falling through a 
Newtonian fluid. We have proven that the equations do not admit oscillations, 
even in the transient, in agreement with general experimental observations. 

From our analysis it appears that the lack of oscillations is due to a delicate 
balance of terms. It is tempting to conclude that an oscillating motion could 
be produced with only a slight modification to the equations. However it is 
important that the solution still remain bounded, and as we have shown there 
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is only one trajectory which is insensitive to the linear instability (3?(a) > 0) of 
the homogeneous equation. Transient oscillations of a falling sphere have been 
successfully modeled by King & Waters using an elastic constitutive model |l5| , 
for which a final steady state velocity is approached. 

In principle, however, one cannot simply modify the differential equation 
( p. 17] ) or even ( 2.15| ) to address the oscillations of a sedimenting sphere in a 
micellar fluid [p. |13|; one must return to the full time-dependent partial dif- 
ferential equation. This was indeed how King & Waters obtained their result 
for a linear viscoelastic constitutive model |[15| , but it is not clear that this ap- 
proach will continue to be fruitful as the complexity of the problem increases. 
Self-assembling wormlikc micellar solutions are thought to have a nonmonotonic 
stress/shear rate relation |22|, po[, based on the existence of an apparently in- 
accessible range of shear rates [|20[ |J. It may be that the dynamics of such a 
nonlinear fluid requires the spatial information inherent in the PDEs, and that 
the ODE reduction discussed here is practically limited to linear models. 
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A Derivation of Equations flZTTTp-flZTT^I) 

The equation describing the transient motion of a falling sphere is 



/ / n 

u + u + \ — 



t „,' 



«'00 



n Jo \i - s 



ds = 1, 



(1.1) 



where u(t) is the velocity of the sphere and k is a non-dimensional parame- 
ter which depends on the relative densities of the sphere and the fluid. This 
integro-differential equation can be converted to a second order ODE through 
the following procedure. If 



(1.2) 



lo \/t-s 

then Abel's theorem (see e.g., (l4|, §3.7]) implies 

4 F(r) 



o Vt-T 



dr = ir [u(t) - u(0)] 



Multiplying (1.1) by 1/ \/t — r, integrating, and using (1.2) yields the equation 



lo Vt=7 



u 

\]t-T 



dr+ I , " dr + ttW - \u(t) - u(0)} = 



1 



dr. (1.3) 
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From (1.1) one observes 

u'(t) 



dr = \ / — (1 — u — u') . 

Vt^T V K 



Substituting this into (1.3) and rewriting yields 
I nt 



1 - 2 



+ («-!)«- ku(0) + 



K Jo Vt~ 



dr. 



(1.4) 



(1.5) 



The desired second order differential equation is now obtained by differentiating 
(|l.5[). In this regard, note that the substitution r = t — x 2 implies 



thus 



J(t) = 



K /" U(t) 



. dr = 2J- 

77 Jo v* — T V 



w(i — x 2 ) dx, 



(1.6) 



Ku(0) 



7T Vt V I" 

Ku(0) 

KU(0) 



ii'(t — x 2 ) dx 



dr 



+ (l-u- u 1 ), 



(1.7) 



where again we have used fll.4f). 



Therefore differentiating (1.5) and using (1/7) yields the second order equa- 



tion 



u" = (k - 2)v! -u + 



l + \l-(u(0)-l) 



(1.8) 



Note that from (1.1) the initial value of u' is determined by the initial value of 
u, i.e., w'(0) = l-u(O). 

Therefore, the equation describing the transient motion of the sphere is 



u" + (2-k)u' + u = 1 + W— (u(O)-l), (1.9) 

V 7rf 

u(0) = f, u'(0) = l-f. (1.10) 
If the sphere starts from rest (i.e., u(0) = 0) then the system reduces to 

u + (2 — k)u + u = 1 — 



t(0) = 0, tt'(0) = 1, 



which is precisely (2.17)-(2.18). 
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